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The connection between quantized vortices and dark solitons in a long and thin, waveguide- like 
trap geometry is explored in the framework of the non-linear Schrodinger equation. Variation of 
the transverse confinement leads from the quasi- ID regime where solitons are stable to 2D (or 3D) 
confinement where soliton stripes are subject to a transverse modulational instability known as the 
"snake instability". We present numerical evidence of a regime of intermediate confinement where 
solitons decay into single, deformed vortices with solitonic properties, also called svortices, rather 
than vortex pairs as associated with the "snake" metaphor. Further relaxing the transverse confine- 
ment leads to production of 2 and then 3 vortices, which correlates perfectly with a Bogoliubov-de 
Gennes stability analysis. The decay of a stationary dark soliton (or, planar node) into a sin- 
gle svortex is predicted to be experimentally observable in a 3D harmonically confined dilute gas 
Bose-Einstein condensate. 



03.75.-b, 05.45.Yv, 42.65.Tg 

Solitons and quantized vortices are fundamental exci- 
tations of non-linear media. Quantized vortices, often 
regarded as an indicator for superfluidity, are topologi- 
cal defects in (2+1)- or (3+l)-dimensional fluids. Dark 
solitons, in their purest form are solitary, nondispersive 
density- notch solutions to (l+l)-dimensional, non-linear 
wave equations with extraordinary stability properties. 
It has been known, however, for many years that solitonic 
wave fronts (also called band solitons or soliton stripes) 
in 2- or 3-dimensional media are unstable Q-g). The 
metaphor of a "snake" instability (SI) has been intro- 
duced in this context by Zakharov H in order to refer to 
the antisymmetric modulation (bending) of the solitonic 
wave- front caused by long- wavelength perturbations [|| . 
Later it has been predicted by numerical studies of the 
time evolution that the SI eventually leads to the forma- 
tion of arrays of vortices with alternating charge [|],[|. 
The first experimental evidence of the SI and subsequent 
formation of vortices was observed in non-linear optics 



More recently, dark solitons have been observed in 
trapped dilute-gas BECs || and the decay of a station- 
ary soliton into closed loops of vortex filaments, much 
resembling smoke rings, has been observed in a spheri- 
cal harmonic trap ||. Stationary dark solitons, like the 
example shown in Fig. 0(a), are nodes (nodal lines or 
planes in 2D or 3D, respectively) in the wavefunction as 
opposed to traveling solitons, which are also referred to 
as gray solitons. Theoretically, the stability of stationary 
solitons in harmonically trapped BECs had been investi- 
gated before by Muryshev et al. [|Io| and Feder et al. [[ll| , 
based on a linear stability analysis using the Bogoliubov- 
de Gennes (BdG) equations. While both papers identify 
a regime of stability for stationary solitons in elongated 
traps at low density as expected from earlier work ||, it 
was conjectured in Ref. [|10[ that the mechanism of insta- 
bility at increasing density was vortex pair production in 
analogy to the SI. Feder et al. O] refined and partially 



corrected the results of Ref. |T(| and predicted the later 
experimentally observed vortex ring formation ||. The 
mechanism of decay at the onset of instability, however, 
has not been fully revealed so far. 




FIG. 1. Two stationary "one-defect" excited states in a 
2D rectangular box trap with hard- wall boundary conditions 
generated by imaginary-time propagation of the NLSE: (a) 
stationary soliton (line node) and (b) stationary svortex state. 
The complex wavefunction ip is represented by the density 
p — | -0| 2 and the phase <j) — arg(^). Each subplot shows a 
surface plot and gray-scale coded plots of the density and the 
phase modulo 2tt. 

In this letter we study the modes of instability of a 
stationary soliton as a function of the transverse confine- 
ment L t , measured in terms of the condensate healing 
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length f US The onset of instability at L t > 6£ is 
initiated by the emergence of a nontrivial stationary state 
[see Fig. 0(b)] of lower energy than the corresponding sta- 
tionary soliton. We call this state a solitonic vortex or 
"svortex". A svortex is a single confined and deformed 
vortex with solitonic properties Q . For transverse con- 
finements of 6£ J$ L t J$ 10£ the strong coupling of the 
stationary soliton to the more stable single svortex is the 
only decay mechanism available in contrast to what has 
been seen and expected in earlier work |^ ]i^]ll| J| . The 
svortex therefore presents the smallest possible unit of 
decay, which persists in geometries where the transverse 
confinement is too tight for vortex ring (in 2D: vortex 
pair) formation. Under less restrictive confinement, 2 and 
then 3 vortex channels open (for L t ^ 10£ and L t }t 13£, 
respectively). The 1,2, and 3 vortex instabilities will be 
seen to correlate perfectly with a BdG stability analysis. 

The essential physics involved reveals itself from study- 
ing the time-dependent Gross-Pitaevskii or nonlinear 
Schrodinger equation (NLSE), which presents the rele- 
vant mean-field theory for a zero-temperature BEC |l4| 
and also applies to non-linear wave propagation in optics 
@: 

id t i/;= [-V 2 + F + ^bH 2 ]^. (1) 

In the dimensionless Eq. ([[]), the condensate wavefunc- 
tion tp(r,t) satisfies the following normalization condi- 
tion: J v ^ |^| 2 <ir= 1, where vb is the volume of a box con- 
taining the trapped condensate. The external trapping 
potential is given by V and g is the nonlinear coupling 
constant. We restrict ourselves to a repulsive nonlinear- 
ity (or defocusing NLSE) g > 0. The relevant size scale 
for nonlinear structures like solitons (l5| and vortices 
is the condensate healing length £ = 1/ ^/gvB\^\ 2 ^ where 
£ = 8iraN/(gvB) is unit of length used in Eq. (|l|) for a 
BEC with N particles and an s-wave scattering length a. 
Note that for fairly uniformly distributed condensates, 
the healing length is given by £ = l/y/g£. The applica- 
tion of NLSE solutions of type of Fig. |l| has been fully 
confirmed empirically by the experiments of Refs. ||[§. 
In tightly confined BECs, the current mean-field theory is 
justified as long as the transverse dimensions are greater 
than £ and £ ^> a is satisfied (see [[l7[). 

We initially consider a 2D rectangular geometry with 
box boundary conditions. The chosen aspect ratio 
length/width of 8 simulates a transversely confined, 
waveguide-like geometry. The stationary vortex-like 
state with a node and phase singularity at the trap center 
was found by imaginary-time propagation and confirmed 
by real-time propagation of the NLSE [[l3|. In addition to 
seeding this relaxation procedure with a suitable phase 
profile, we also restricted the symmetry of the density 
\ijj\ 2 to be even in both spacial directions. A second sta- 
tionary state (a dark band soliton) was also generated by 
imaginary-time propagation with the constraint of odd 
symmetry in the longitudinal direction of the trap. Fig- 
ure 111 shows the resulting wavefunctions. The vortex-like 



wavefunction of Fig. |T|(b) is clearly distorted and affected 
by the tight traverse confinement of 8£. We have argued 
in Ref. |lj that such a tightly confined vortex acquires 
solitonic properties and therefore should be called soli- 
tonic vortex or svortex, further discussed in [[18|. 
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FIG. 2. Svortex and soliton properties as a function of the 
transverse confinement Lt. Part (a) shows the excitation en- 
ergies of the stationary soliton and svortex state. The energy 
difference is shown on a logarithmic scale in part (b). The 
simulations were done in a rectangular box of size 8£ x 64£ 
for different values of the nonlinear coupling constant g, which 
changes the effective transverse confinement L t /£ = Sy/g. 

Figure ||(a) shows the excitation energy of the station- 
ary svortex and soliton state with respect to the ground 
state as a function of the transverse confinement L t in 
terms of £. For a given wavefunction -0, the energy is 
given by the formula E = /^ b (-^*V 2 ^ + #v B /2|^| 4 + 
V^l 2 ) dr. The soliton excitation energy exhibits a linear 
dependence on the variation of the length scale reflecting 
the localization of the excitation in one and extension in 
the other spacial dimension. The svortex energy is always 
lower than the soliton energy and grows more slowly with 
the box size, reflecting the expected logarithmic behavior 
for large boxes p6| . Below a critical confinement corre- 
sponding to a box width of w 6£ we do not find any sta- 
tionary svortex solutions but instead the imaginary-time 
propagation converges to the soliton solution. The loga- 
rithmic plot of the energy difference between the soliton 
and the svortex energies shown in Fig. ||(b) very much 
indicates a nonanalytic curve joining or curve crossing. 
Following the svortex solution from wide confinement to 
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the critical point, the vortex wavefunction shows an in- 
creasingly deformed density and squeezed phase signa- 
ture (see Fig. |l|) and eventually coincides with the soliton 
wavefunction at the critical confinement. 

The band solitons and svortex states from Fig. ^|(a) are 
stationary states. In wide enough confinement, however, 
band solitons may exhibit the SI mentioned earlier: Tiny 
imperfections of stationary band solitons may lead to a 
transverse modulation and grow during real-time propa- 
gation at an initially exponential rate. The stability of 
stationary solutions of the NLSE can be tested in a linear 
stability analysis employing the famous BdG equations 
]l9| , which can be derived from a linear-response expan- 
sion of the time-dependent NLSE pQ[ . In the units of 
Eq. |l| the BdG equations read: 

Cuj{r) - gv B [^(r)] 2 Vj(r) = ejU^r)j, (2) 
Cvj(r) - gv B [ip*(r)] 2 Uj (r) = -ejVj(r), (3) 

with C = -V 2 + V(r) + 2gv B \ifj(r)\ 2 - /i, and \i is 
the chemical potential of the stationary wavefunction 
ip(r, t) = exp(— ifittyfr). The solutions of the BdG equa- 
tion with eigenvalues ej and eigenvectors (uj,Vj) have the 
following interpretation in terms of small-amplitude mo- 
tion around a stationary solution of the NLSE pM : Small 
positive €j at positive "norm" rjj = J(\uj\ 2 — \vjf)dr de- 
scribe small oscillations around the stationary state with 
increasing energy. Solutions with negative eigenvalue ej 
and positive rjj are called anomalous modes. They indi- 
cate a continuous transformation of the stationary state 
to a state of lower energy. Anomalous modes exist for 
the trapped vortex as well as for dark solitons in ID and 
merely express the thermodynamic instability of these 
excitations. Complex or purely imaginary eigenvalues ej, 
however, indicate a dynamical instability. They further 
imply rjj = 0. 

Figure || shows the purely imaginary and anomalous 
eigenvalues of the BdG equation for a stationary_band 
soliton in a rectangular box of dimension b x 16£ as a 
function of the box width b « L t at constant density. For 
narrow traps with b ;$ 5.5£, we find one anomalous but 
no complex eigenvalues, like for ID solitons. Addition- 
ally, the soliton wavefunction shows no appreciable decay 
in real-time propagation seeded with noise (see insets). 
Also collisions of noisy gray solitons show the robust, 
particle-like behavior expected from ID soliton theory 
[p~8|| . For trap widths 5.5£ J$ b ^ 9.5£ one purely imag- 
inary eigenvalue exists in the BdG eigenvalue spectrum. 
According to the numerical results, the emergence of this 
imaginary eigenvalue coincides with the emergence of the 
svortex as a symmetry-breaking stationary state of lower 
energy than the corresponding band soliton. Increasing 
the box width, a second and eventually a third imaginary 
eigenvalue appears. The stability of the stationary soli- 
ton was probed using real-time propagation seeded with 
a small amount of white noise |2l|] . While there is no 
appreciable decay in tight confinement, we clearly find 
that the soliton instability is associated with the forma- 



tion of one, two, and three vortices in the regimes where 
one, two, and three imaginary eigenvalues are present as 
shown in the insets of Fig. ||. The BdG eigenvectors Uj, 
localized within about one healing length from the nodal 
line of the soliton, also support this result . The pat- 
terns shown in Fig. || are by no means stationary but 
rather form transient states followed by incomplete re- 
currences of the nodal line and eventual further decay 
where vortices move to edge of the trap and vorticity is 
destroyed. The complicated dynamical patterns showing 
a mixture of decay and strong mode coupling are cer- 
tainly due to energy conservation in the NLSE and to 
the small scale of the trap used in the simulation where 
radiated phonons linger. We expect further stabilization 
of the vortex patterns in longer traps where energy re- 
leased in the decay process can distribute itself over a 
larger area. The observed decay patterns vary depend- 
ing on the exact form of the initial perturbation by noise. 
In contrast to the soliton, the stationary svortex shows 
an entirely real BdG spectrum with one anomalous mode 
also shown in Fig. |[ Further, real-time propagation of 
perturbed svortices shows no appreciable decay. In this 
sense, the svortex is the more stable object than the sta- 
tionary soliton. 




width b/^ 

FIG. 3. Eigenvalues in the BdG spectrum of the station- 
ary soliton in a 2D box as a function of the box width at 
constant average density. The insets relating to box widths 
of b = 7£, 12£, and 16£ show density plots of transient patters 
in the decay of the perturbed soliton state [21] after real-time 
propagation for t — 26, 31, and 26 in the units of Eq. |l|, 
respectively. The perturbed soliton at b — 5£, on the con- 
trary, shows no appreciable decay after 100 time units. The 
imaginary modes are marked according to the nature of the 
eigenvector u leading to single- vortex (A), double- vortex (o) 
or triple- vortex (v) decay. The anomalous mode of the nodal 
plane-state (x) and stationary single- svortex state i^k) are 
also indicated. 

Finally, we would like to comment on the 3D harmonic 
trap geometry studied earlier in the JILA experiment |J 
and in theoretical work by Muryshev et al. |l(J and Feder 
et al. |ll[. Both experiment Q and theory |ll|] report 
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vortex-ring formation during the decay of a stationary 
soliton (nodal-plane state) in spherical fllTTf and elon- 
gated jn| geometries at fairly high densities, which is 
indicated by the nature of complex modes in the BdG 
spectrum of the stationary soliton jnj]. It has also been 
pointed out that the BdG spectrum becomes entirely real 
at sufficiently low particle number or high aspect ratio in 
elongated traps. However, the decay mechanism in the 
presence of a single imaginary mode in the BdG spec- 
trum (as shown in Fig. 4 of Ref. |ll|]) of a stationary 
soliton in an elongated trap has not been revealed so far. 
Imaginary- and real-time propagation clearly show that 
a stationary svortex solution exists in this regime and 
that it has lower energy than the stationary soliton. The 
density and phase profile of the stationary svortex state 
are very similar to the dynamically generated pattern 
shown in Fig. |^. This figure shows the transient decay 
product of a perturbed stationary black soliton after real- 
time propagation for 100 ms, seeded initially with 0.01% 
white noise. The parameters of this 3D harmonic trap 
(N = 10 4 atoms of Na) correspond to Fig. 4 of Ref. |j]] 
at an aspect ratio of oo p /co x = 4 with uo x = 2tt • 50rad/s. 
The corresponding imaginary Bogoliubov mode Uj has an 
azimuthal coordinate dependence of exp(z0), where <\> is 
the azimuthal angle, and much resembles the first imag- 
inary BdG mode in the 2D box discussed above. The 
predicted decay of the band soliton into a single svortex 
has not been seen, or predicted, before and should be 
easily observable with current experimental techniques. 




FIG. 4(color). Svortex in a 3D elongated harmonic trap 
generated by decay of a perturbed stationary soliton state. 
Shown is the surface of constant density (at 0.16 of the max- 
imum density) and a color-coded plot of the phase in the 
horizontal plane intersecting the trap center. The transverse 
confinement Lt/£ ~ 7.7 was computed by the maximum value 
of the line integral J c £(r) -1 ds taken along the transverse di- 
mension, which is more appropriate for measuring the trans- 
verse confinement of inhomogeneous condensates than the box 
width [13]. For details of the simulation see text. 

Concluding, we have identified the fundamental modes 



of the SI for transversely confined geometries: Produc- 
tion of 1, 2, and 3 vortices correlates with imaginary 
modes in the BdG eigenvalue spectrum. Departure from 
the quasi- ID regime of stability of solitons is indicated 
not only by a linear stability analysis but also by the 
emergence of a solitonic vortex or svortex as a stationary 
state of lower energy than the corresponding dark soli- 
ton. We demonstrated that the decay of a soliton into a 
single svortex is a fundamental mode of instability in 2D 
box geometry and 3D elongated harmonic traps. 
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